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Appendix A 

Three aspects of the invention are: 

1. The recognition of the desirability of displaying to a financial investment customer in real 
time, for example on a World Wide Web site, the probability distribution governing the 
price of a particular asset (e.g., a stock) at a selected future time. 

2. The recognition that such probability distributions can be derived from option prices for 
that asset, or for related assets, which are readily available in real time. 

3. The recognition that probability distributions involving several asset prices simultaneously 
are useful to investment customers in several contexts, especially in exploring hypothetical 
scenarios, and that single asset distributions such as (but not restricted to) the above can 
be meaningfully incorporated into multivariate distributions, manageably determined. 

In this appendix we first describe a basic method for deriving probability distributions for 
single assets from option prices. We next describe improvements on this basic method to address 
various practical issues. Then we take up the multivariate case and show how to extend this 
kind of single asset price distribution, or any other, to the multivariate case. Finally, we consider 
a number of novel multivariate applications, with emphasis on scenario exploration. 

1 Basic method 

A call option is an option to buy an asset (e.g., a stock) at a certain price x (called the strike price) 
on a given expiration date T days in the future. (An option exercisable only on the expiration 
date is called a European-style option; for simplicity we will consider in this discussion only this 
type of option. 1 ) Similarly, a put option is an option to sell an asset at a strike price x on a 
given expiration date. (The "European-style" assumption of no possible early exercise is more 
important here, but can also be ignored for puts that axe not too deeply "in the money." ) 

Let c(x) denote the price of a call option on an asset at strike price x, and p(x) the price of a 
put option. Such prices are established by options market-makers. We have realized that such 
prices implicitly contain information about a "market view" of the probability distribution of 
the price of that asset at the expiration date. 

In a simple but precise form, this market view can be stated as follows. Suppose that we were 
given the call price curve c(x) or the put price curve p(x) as a continuous function of the strike 
price x for all x > 0. Then, the second derivative of either the call or the put price curve is the 
market view of the risk-neutral probability density function (pdf) f(x) of the asset price at the 
expiration date. In other words, f(x) = d f (x) =p ,f (x). 

The idea that option prices determine some kind of implied probability distribution is fairly 
well known in the financial literature. The idea that a pdf can be computed by taking the 
second derivative of a continuous option price curve is known in the academic literature, but it 
does not appear to be very well known. For example, the standard textbook "Options, Futures, 

*Even allowing for possible early exercise, most liquidly traded call options without large dividends can be 
treated as if there were no possibility of such exercise, since sale of the option is usually a better alternative; 
therefore, these call options behave similarly to European-style options. 






and Other Derivatives," by John C. Hull (Fourth Edition, 1999; Prentice-Hall) mentions implied 
probabilities, but not the second-derivative method. The best reference that we have been able 
to find is J. C. Jackwerth and M. Rubinstein, "Recovering probability distributions from option 
prices," J. Finance, vol 51, pp. 1611-1631 (1996), which has only six prior references. 

The risk-neutral distribution (at a fixed future time T, for a fixed asset) is defined as the price 
distribution that would hold if market participants were neutral to risk, which they generally 
are not. However, many asset pricing theories, such as those underlying Black-Scholes option 
theory and most of the variations found in the Hull book above, allow for the true risk-averse 
asset price distribution to be obtained from the risk-neutral distribution /(x) just by adjusting 
the latter by an appropriate risk premium: If there are no dividends, the true distribution is just 
/(xe^" r ) r )), where fi - r is the expected annual return rate for the stock in excess of the risk 
free rate r. We use a variation on this simple format, slightly modified to allow for dividends (see 
below), though our invention could also work well with a more complicated adjustment. In this 
format, a value for ji — r must still be supplied. We use as a default the "consensus estimate" 
taken from the textbook "Active Portfolio Management" (1995) by Grinold and Kahn. These 
authors note a long-term average value of the risk premium to be 6% per year, and suggest 
multiplying this number by the stock's beta to get fi — r. The parameter beta is the slope of the 
line giving a regression of the stock in question against a market portfolio, often taken as the 
S&P 500. This is the well-known CAPM estimate for the expected excess return. Whether good 
or bad, its stature as a consensus estimate makes it suited to our aim of providing a market 
view, though it is only a default. Our invention, which provides the risk-neutral component 
of the probabilities, could work with other estimates for the risk-averse adjustment parameter 
\i — r and with any explicit scheme for adjusting the risk neutral probability density to the 
risk-averse probability density. It is worth pointing out that, for shorter time periods-even a 
month or two-the risk adjustment required is small and generally overwhelmed by fluctuations 
in the risk-neutral distribution itself. 

We give a brief proof that the second derivative procedure gives the correct risk-neutral prob- 
ability distribution. As in Hull, we may calculate the European call or put price as an expected 
value in the risk-neutral distribution. 

If the actual value of the asset on the expiration date is v, then the value of a call option at 
strike price x is max{x; — x,0}, and the value of a put option is max{x — v,0}. If the actual 
value is a random variable with pdf f(v), then the expected value of a call option at x at the 
expiration date is 



The current values c(x) and p(x) may be obtained by discounting ct(x) and pr{x) by e~ rT , 
where r is the risk-free interest rate, but for our purposes, forecasting probability distributions 
at time T, we do no discounting, and henceforth just write c(x) = cr{x), p(x) = pr(x). 

Parenthetically, from these expressions we observe that 



c T (x) 




and the expected value of a put option at x at the expiration date is 



p T (x) = E v [max{x - 0}] = / (x — v)f(v) dv. 

Jo 




p(x) — c(x) = 




expected value of the asset at the expiration date 



where 5* = E v [v] is the expected value of the asset at the expiration date under the risk-neutral 
distribution. (If there are no dividends, then 5* = se rT ; if there are dividends, then in general 
it is necessary to subtract from se rT the value at time T of the dividends.) This well-known 
relation is called put-call parity; it shows why either price curve carries the same information. 

Prom the above expression for c(x), it follows that its first derivative is 



c f (x) = - / f(x) dx = F(x) - 1, 

Jx 



where F(x) = f£ f(v) dx is the cumulative distribution function (cdf) of the random variable 
v. To prove this, note that v — x = f£ dz. Therefore 

poo roo rv poo poo poo 

c(x)= (v - x)f(v) dv = / dv I dzf{v)= dz dv /(«) =/ dz(l-F(z)), 

Jx Jx Jx Jx Jz Jx 

where we interchange the variables v, z to integrate over the two-dimensional region 72. = {(v y z) : 
x < z < v}. The last expression implies that d(x) = —(1 — F(x)). 

Prom put-call parity, it follows similarly that 

p'(x) = l+c'(x)=F(x). 

Since the cdf and pdf are related by F*(x) = /(#), these expressions in turn imply that the 
second derivative of either c(x) or p(x) is the pdf f(x): 

c"(x)=p"(x) = F'(x) = f(x). 
The general character of the option price curves c(x) and p(x) is therefore as follows: 

• For all x less than the minimum possible value of v (i.e., such that F(x) = 0), c(x) = 
E v [v] — x = s* — x and p(x) = 0. In other words, c(x) is a straight line of slope —1 starting 
at c(0) = E v [v] ~ 5*, while p(x) = 0. 

• For all x greater than the maximum possible value of v (i.e., such that F(x) = 1), c(x) = 0 
and p(x) = x — s*. In other words, p(x) is a straight line of slope +1 and x-intercept 5*, 
while c(x) = 0. 

• These two line segments are joined by a continuous convex U curve whose slope increases 
from —1 to 0 for c(x), and from 0 to +1 for p(x). 



We note that the fact that the mean E v [v] of the pdf f(x) is 5*, the value in future dollars at 
time T of the underlying price s (less the value of any dividends), implies that option prices must 
be constantly adjusted to reflect changes in the underlying price 5, even if there is no market 
activity in the options. 

The fact that s* = E v [v] also implies that an option price curve can make no prediction about 
the general direction of the underlying price s. However, the option price curve does predict the 
shape of the pdf /(x), and in particular its volatility. 
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^ased cm a finite subset of bid-ash^f c 



1.1 ^prcDdxratioais cased cm a finite subset of bid-asked option prices 

In practice, option prices c(x) and p(x) axe quoted only for a finite subset of equally-spaced 
strike prices x, namely x = nA for integer n and spacing A. We denote c(nA) and p(nA) by On 
and p n , respectively. Moreover, quotes specify only a bid-asked spread, not exact prices. In this 
subsection we give methods for dealing with these problems. (Most of the Jackwerth-Rubinstein 
paper (op. cit.) is concerned with these kinds of curve-fitting problems.) 

The first derivatives d{x) and p f (x) at x = (n + |)A may be estimated by the first differences 

% _ Cn+l ~Cn - _ Pn+1 ~ Pn 

The corresponding estimates of the cdf F n+ i = F((n+ j)A) are 



1 . 



Thus, using both bid and ask prices for both European-style puts and calls, one can compute 
four different estimates for the cdf F n+ i , which can then be combined into a single estimate. 

This combination will preferably take into account whether x = (n + ^-)A is much less than 
the underlying price 5 ("deep out-of-the-money"), near s ("near the money"), or much greater 
than s ("deep in-the-money"), according to the different patterns of setting bid-asked spreads in 
these different ranges. Another consideration is avoiding quotes near prices where early exercise 
is likely, such as deep in-the-money puts. 

Similarly, the second derivatives d f (x) and p n (x) at x = nA may be estimated by the first 
differences of the estimates of the first derivatives; e.g., 



d^ 1 - d 

c n — 



n+| _ Cn+l ~ 2Cn + Og-i 



A A 2 

We may take d f n or p n n or some combination as above as our estimate f n of the pdf /(nA). 

Note that since f(x) > 0, option prices should satisfy a convexity condition, e.g., Cn+i — 2cn + 
Cn-i > 0 for call option prices. Indeed, violation of this condition would allow making money 
via a risk- free "butterfly straddle" involving buying one call option at (n + 1)A and another at 
(n — 1)A, and selling two call options at nA. A similar result holds for put options. 



1.2 Dynamic estimates 

The methods considered in the previous subsection allow estimation of the cdf and pdf at a 
subset of A-spaced values of x, based on a static set of option quotes at a particular time. 

As previously noted, however, option prices must change continually in response to changes 
in the underlying price s. Let s* denote the corresponding forward price at expiration (the 
price 5 evaluated with interest). Suppose this price (measured in dollars at expiration) moves 
up (or down) by a small amount, an increment e in its logarithm, say, with little or no change 
in volatility. Here e may be viewed as, approximately, the percentage move 6/s* caused by a 
move of 6 in the (forward) stock price. We expect in this situation that (forward) probability 
distribution for the stock price will just be shifted by e in the log domain. That is, the distribution 
will appear to be identical there, except with a mean shifted by e. Thus, the value of the new cdf 
at x = e* nx is F(e^ nx ~ e ^) — F(x/a) y where F denotes the original cdf with distribution mean 
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$*, and a = e 6 . A reasonffole call price functional equation that gives the same effect, upon 
differentiation, is 

ac(s*)X/a) — c(as*,x), 

where c(s*,x) denotes the price, in dollars at expiration, for a call option at strike x when the 
underlying is at price s. Note in this equation that all other variables, such as volatility, are 
assumed to be the same, which will only be approximately true, even for very small values of €. 

But, assuming this approximation, we can think of an option price at strike x ) measured when 
the (forward) price has moved to as*, for a near 1, as giving instead a times the price of an 
option at strike x/a, but corresponding to the current underlying price s*. Considering all 
the strikes at which options axe frequently quoted, and thinking additively, we can effectively 
observe c(x) (and p{x)) for a different subset of approximately equally-spaced strike prices, 
roughly x = nA — S for various values of 5 = e$*. Some care must, of course, be taken to ensure 
simultaneity of prices, of option and underlying. For this reason, we may prefer to consider the 
values of nA (corresponding to the various standard strike values) separately, and synchronize 
observed time of sales for an option at a given strike with the underlyling security. Implied 
volatilities (discussed below) could be monitored, to ensure their changes relative to e were 
small. 

Using a similar technique to that described in the above paragraphs, meaningful average option 
prices for a given strike can also be computed, using thin strike intervals and using either short 
time intervals or time series methods (time averages weighting the present more than the past). 
Note that, without the framework described in this subsection, the computation of "average" 
option prices at a given strike are problematic when the stock price varies in the period over 
which the average is taken. 

To summarize: Given enough movements of the underlying price, we can effectively observe 
prices and compute estimates as above for a much more finely quantized subset of strike prices 
x, and provide a framework for improving accuracy through averaging methods. 



2 Methods for extrapolation and smoothing 

There are two major limitations to the basic methods of the previous section. One is that option 
quotes are available only for certain expiration dates. Another, not so obvious, is that option 
quotes are reliable primarily for options in which there is substantial market activity. These 
would typically be nearer- term options at strike prices near the money (the underlying price). 

To extend our prediction methods to times other than expiration dates and over wider ranges 
of strike prices (and also to help reduce "noise" in our displays), we use extrapolation and 
smoothing techniques. We have found that it is advantageous to do extrapolation and smoothing 
in the volatility domain. 

There are many reasons for this advantage. For example, option practitioners axe well aware 
of the kinds of shapes that the volatility curves (sometimes called "volatility smiles" ) have had 
historically, in various markets, and how these curves vary with time; this can be a guide to 
imposing structure on the smoothing curves to prevent overfitting of possible artifacts. Many 
records have been kept of the volatilities implied by option prices, and it is easy to examine how 
in the past they have changed with respect to price behavior. For example, the Chicago Board of 
Options Exchange makes public its average near-the- money volatility index (now called vix) for 
S&P 100 options back to 1986. Finally, it is easier to work visually with volatility curves, which 



would theoretically be flat^^f (x) were lognormal, than with visual differences in near-lognormal 
pdfs, which can all look very much alike. Mathematically, model improvements can be made 
in the volatility domain just by changing coefficients of low-degree polynomial approximations, 
even though these affect higher-order terms in power series for the corresponding cdfs or pdfs. 

The following subsections explain more precisely how to work in the volatility domain. 



2.1 Lognormal pdfs 

The standard Black-Scholes theory of option pricing (see Hull, op. ext.) yields a lognormal pdf 
f(v) whose expected value is E v [v] = $*, such that In v is a Gaussian (normal) random variable 
with variance <r 2 T, where the parameter a is called the volatility rate of the asset, and T is 
the time to expiration. By a standard property of lognormal distributions, this implies that the 
mean of \nv is E v [ln?;] = Ins* — cr 2 T/2. 

Prom this pdf follows the famous Black-Scholes call option pricing formula [Hull, Appx. 11 A]: 

c(x) = E,,[max{v - x,0}] = s*N(di(x)) - xN(d 2 (x)), 

where N(d\ (x)) and N(d 2 (x)) are values of the cumulative distribution function of a Gaussian 
random variable of mean zero and variance 1 at the points 

d 2 (x) +<jvf ; 

E v [lnt>] — lnx 
ay/f 

(Recall that our version of the call price is not discounted, and is given in dollars at time T, and 
that 5* is today's stock price, valued in dollars at time T, less the value of any dividends.) Note 
that cry/T is the standard deviation of lnt>; therefore — d 2 (x) is just lnx, measured in standard 
deviations from the mean E v [lnv]. 

Similarly, by put-call parity, we have the Black-Scholes put option pricing formula 

p(x) = c(x)+x- s* = s*{N(d x (x)) - 1) - x{N(d 2 (x)) - 1) = xN(-d 2 (x)) - s*N(-d 1 (x)). 



d\ (x) = 
d 2 (x) = 



ln(s*/:r) + q 2 r/2 

ln(s*/x)-a 2 T/2 
ay/T 



Taking the derivative with respect to x, and using s*N*(d\ (x)) = xN'{d 2 (x)) and d[ (x) = d l 2 (x) 
(the latter equation holding under the assumption of constant volatility, which we will later 
drop), we obtain 

F(x) = c'(x) + 1 = -N(d 2 (x)) + 1 = N(-<k(x)). 

Now F(x) is the probability that v < x, which is equal to the probability that Inv < lnx, which 
since Inv is Gaussian with mean E v [lnv] and standard deviation a\ff is given by 

F{x) = Pr{v <x} = Pr{lnv < lnx} - N ^ nx -^ lnv } ^ = ;v(-d 2 (x)). 

Thus we have verified that the Black-Scholes pricing formulas give the correct cdf F(x). The 
derivative of F(x) will thus yield the correct lognormal distribution /(x) = F f (x). 





2.2 Giaracterizatidr^jf general odfe 

Now let F(x) be an arbitrary cdf on R+; i.e., a function that monotonically increases from 0 to 
1 as x goes from 0 to infinity. For simplicity we will assume that F(x) is strictly monotonically 
increasing; i.e., f(x) = F'(x) > 0 everywhere. Then there exists a continuous one-to-one 
"warping function" y : R+ — ► R such that F(x) = N(y(x)) everywhere; i.e., such that the 
probability that a random variable v with cdf F(x) will satisfy v < x is equal to the probability 
that a standard Gaussian random variable n with mean zero and variance 1 will satisfy n < y(x). 
Similarly, there is an inverse warping function x(y) such that F(x(y)) = N(y). 

Given the warping function the cdf F(x) may be retrieved from the relation F(x) = 

N(y(x)). Therefore the cdf F(x) completely specifies the warping function and vice versa] 
i.e., both curves carry the same information. 

If F(x) is the cdf of a lognormal variable v such that Inv has mean E v [lni;] = Ins* — a 2 T/2 
and variance a \ — a 2 T, as in the previous subsection, then the warping function is given by 



For this reason we may sometimes write y(x) as — ^(x), even when the cdf is not lognormal so 
that the right-hand equation above for d 2 {x) does not hold. 

2.3 Implied volatilities 

If f(x) is not lognormal, then the Black-Scholes pricing formulas do not hold. Nonetheless, given 
an option price c(x) or p(x), it is common practice to define the implied volatility <r(x) as the 
value of a such that the Black-Scholes pricing formula holds, for a given x, s and T. 

The implied volatility curve a(x) so defined is a function of the strike price x, which is constant 
if and only if the pdf f(x) is actually lognormal. In practice, it is typically a convex U curve, 
called a "volatility smile." See, e.g., Hull, Chapter 17. 

From Subsection 2.1, we can see that there is a second method of calculating implied volatilities, 
as follows. Suppose that we have an estimate of the cdf F(x). Define the cdf-implied volatility 
cri(x) as the value of a such that the Black-Scholes cdf formula F(x) = N(—d2(x,a,T)) holds, 
for a given x, s and T. 

The first method has the advantages of being defined directly from raw price data, and of being 
well understood in the financial community. However, the second method has the following 
advantages: 

1. It is easier to calculate, at least from estimates of F(x)\ 

2. It gives a simpler and arguably more intuitive relationship between volatility and the cdf 
F(x). If we use the traditional implied volatility a(x), then the relationship is instead 



y (x) = -d 2 (x) = 



lns-(lns*-q 2 r/2) _ lnx - E v [lnv] 



F(x) = N(-d 2 (x)) + —a'(x). 



3. It fits better with the multivariate theory to be developed below. 



We have observed that the Two curves a(x) and a\ (x) seem to be fai^y similar, at least as to 
the direction of their slope, and axe generally not too far apart in value "near the money" . Also 
cti(x) = <r(x) whenever a(x) has zero slope, though <j\{x) is a little smaller than a(x) when the 
slope a'(x) is negative (which often occurs for stocks). See the above equation. Finally, one 
function is as ad hoc as the other. Therefore, because of the above reasons, we generally prefer 
to use the cdf-implied volatility curve <ti(x). 

In any case, it is clear that either a(x) or <j\ (x) contains the same information as any of the 
curves c(x), p(x), F(x) or f(x). Prom cr(x) or cri(x) we can recover c(x) or F(x) using the 
Black-Scholes call option pricing or cdf formula, and from this we can obtain all other curves. 



2.4 Extrapolation and smoothing in the volatility domain 

The volatility curve a(x) or <j\ (x) may be calculated pointwise from the corresponding curve 
c(x) or p(x) to give a set of values at a finite subset of strike prices x. Each of these values may 
be deemed to have a certain degree of reliability. 

It is then a standard problem to fit a smoothed and extrapolated curve a(x) or oi(x) to these 
points, taking into account their relative reliabilities. Any standard smoothing and extrapolation 
method may be used. In general, the usual problems of avoiding overfitting or oversmoothing 
must be addressed. 

It is well-known that implied volatilities also vary with time. We generally wish to estimate 
curves cr(x,T) or <7i(x,r) as replacements for the constant volatility a in the Black-Scholes 
formulas, e.g., c(x) = c(x,cr,T) or F(x) = N(-d2(x } a,T)). 

In an especially meaningful example, we have experimented with a class of smoothing algo- 
rithms used in "Implied volatility functions: Empirical tests," by B. Dumas, J. Fleming and 
R. E. Whaley, J. Finance, vol. 53, pp. 2059-2106, Dec. 1998. These authors fit an implied 
volatility curve <r(x), for the purpose of setting up a "strawman" option price model for testing 
(and defeating) a theory regarding the role of volatility in option pricing. Their "strawman" 
option pricing model c(x) was obtained by putting the resulting smoothed curve back into the 
Black-Scholes call formula. It is a "strawman" ad hoc model, because no intuitive notion of 
stock volatility could possible vary with strike price, which the stock never "sees." Nevertheless, 
their model performed admirably, surpassing in predictive power the highly regarded "implied 
tree" method. One possible explanation offered was that their model mimicked in a smooth way 
interpolation methods actually employed by practitioners in the options markets. (See the dis- 
cussion of "Volatility matrices" in Hull, cited above.) Such an approach to option pricing seems 
ideal to us, because of its accuracy and because its underlying rationale represents a market 
view. Thus, we use the Dumas-Fleming- Whaley model for our own entirely different purpose, 
that of forecasting probability distributions. All that is necessary is to differentiate their call 
price model, which, conveniently for us, is a smooth function of strike price and other standard 
variables such as time, current stock price, and the risk- free rate of interest. The formula for 
the cdf F(x) is, as before, this derivative with 1 subtracted, or 

F(x) = N(-d 2 (x)) + ^o'(x). 

We can make this very explicit. We have 

dc 

7 r- = xVTN'(-d 2 (x)) 
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where N f (z) denotes standard normal density, while cr'(x) may be computed by differentiating 
the Dumas-Fleming- Whaley fitted volatility curve. The latter has the form 

a(x,T) = ao + a\x + a 2 x 2 + a^T + a 4 T 2 + a 5 xT. 

The coefficients {a*} are determined by regression. This kind of quadratic curve-fitting is easily 
implemented. Dumas-Fleming- Whaley impose a constraint to prevent their volatilities from 
going below 0 (or even below 0.01), and we have imposed further constraints on extrapolations 
(which we often carry out beyond the range of their tests), to ensure that the final cdf does not go 
below zero or above one. We have experimented with other variations on their basic approach, for 
example, using linear interpolation in the time domain, where we do not need to take deivatives. 
Our methods would, of course work, with any approach, possibly quite different, to volatility 
curve-fitting, though the general Dumas-Fleming- Whaley approach has many things going for 
it: accuracy, conformity to marketplace use of Black-Scholes, smoothness (differentiability, in 
particular), conformity to historical experience regarding the smile structure of volatility curves 
(especially important for extrapolation), and simplicity (which, beyond ease of implementation, 
helps avoid overfitting). These advantages are achieved in a probability context that was not 
considered in the paper where these volatility curves were introduced. 



3 The multivariate case 



fft The methods in the previous sections are capable of generating a display of raw or smoothed 

^ and extrapolated probability distributions for any optionable asset. Option prices are quoted 

U on a large number of securities, as well as on certain indices, such as the S&P 500. 

^ However, an investor would also like to know future probability distributions for: 

£3 

• His or her entire portfolio; 

C3 • Mutual funds; 

£3 

is, • A security without a quoted option; 



£3 



A security in a hypothetical scenario. 



All of these questions involve considerations of several securities at once, and the probabilities 
of their simultaneous configuration of prices. This is clearly a consideration in the first two 
items above, but also enters in the third, where we would want to extract as much information 
as possible about the security without a quoted option price from those correlated with it that 
do have quoted options. Finally, in scenario analysis there are many questions that involve 
considering the probabilities of several security prices occurring at once, including changes in 
factors influencing the market that might be modeled by changes in a portfolio of those securities 
most affected. We will take all of these issues up in the remainder of this document, but for now 
we just try to give a basic introduction. 

For a portfolio of securities, or a mutual fund, we are interested in a composite asset of the 
form 

x = h\x\ + h 2 x 2 + ■••+■ h n x n , 

where the x\ are all assets for which we individually know the cdf F(xi) or the pdf f(xi). To give 
our method the most flexibility, we do not require that this knowledge come from any particular 
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procedure, though we favo^me approach of the preceding two sections^Iowever, even for some 
securities or indices with a quoted option, we might not feel there was sufficient option activity 
to justify a full fit of a volatility curve, and might take a cruder substitute, even a flat straight 
line based on an average of available implied volatilities. In addition, it is convenient to allow 
the possibility that a few assets we axe monitoring might not have any quoted option at all; this 
is easily accommodated by, say, using a flat volatility curve with a historical value for volatility. 
For testing purposes and comparisons we might even want to consider a list of assets with all 
volatility curves given this way. In any case our methodology here is very general, and we only 
require that we know warping functions yi(xi) such that F(xi) — N(yi(xi)) for all i. If the 
asset has an active options market, then the warping functions may be determined by either 
first estimating F(xi) directly from (finite differences of) options price data, as in subsection 
2.2, or by using the approach discussed later in Section 2 of extrapolating and smoothing in the 
volatility domain. In the latter case we have an explicit form of the warping function yi(xi) in 
terms of a fitted volatility curve a\(xi,T) as yt{xi) — — d2(xi,ai(xi,T),r), and this equation 
can also be used with any volatility curve with the assets above, that might have fewer or no 
traded options. In a later section we will discuss portfolios in the logarithm domain, possibly 
containing long and short positions. One can think of warping them to standard normal directly, 
subtracting the mean and dividing by the standard deviation. Alternately, to keep our notation 
uniform, one can invent an asset with price xi such that — ^(^i) gives this warped value (using 
for <T\{x) the observed historical volatility). But we wish to emphasize that the method we are 
describing works with ANY single- variable warping functions, even using a different one for each 
variable. The only further substantitive ingredient is the plausibility of using JOINTLY normal 
distribuitions, which we now discuss. 

The general problem is to find a multivariate probability distribution for the complete set of 
variables (xi, . . . , x n ), or equivalently for their logarithms. In simple financial models generaliz- 
ing the Black-Scholes framework, the multivariate distribution of the logvariables is multivariate 
(i.e. j jointly) normal; see Musiela and Rutkowski's book "Martingale methods in financial mar- 
kets" (1999). This implies that all portfolios of these logvariables are jointly normal, and can 
also be used with other logvariables and portfolios of them to form a jointly normal distribution. 
Thus, if we wish, it is reasonable to use BARRA (or functionally equivalent) factors as single 
(log)variables in our model, using, say, individual normal distributions for them based on histor- 
ical volatility. These factors may represent fundamentals of companies or even macroeconomic 
variables such as interest rates. We do not further discuss such factors, but refer to the book 
of Grinold and Kahn cited above, which also describes how to closely approximate them as 
portfolios of security returns. Our preference is to not use BARRA factors directly, but stay as 
much as possible in the world of optionable securities, and address questions involving BARRA 
factors in terms of approximating portfolios consisting mostly of optionable securities. (But for 
testing and comparisons, it is still useful to be able to include them directly, and we do have 
that capability.) 

Now we certainly do not wish to use only the simple multidimensional Black-Scholes model, 
which would not directly allow the nonlognormal input from our single-variable distributions 
based on the options markets. At the same time, option prices on individual assets do not tell 
us anything about how assets interact, in particular, their correlations. Fortunately, correla- 
tions may be estimated from past (historical) data, and may be viewed as covariances for data 
that has been standardized (has standard deviation 1). Each multivariate normal distribution 
is determined by its mean and covariance matrix. Thus, a natural approach is to use the in- 
dividual distributions to transform or "warp" the variables to standard normal, then impose a 
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multivariate normal stru^^ based on the correlation matrix. ThiFprocedure is independent 
of the individual warping functions, which may be different for different individual variables, 
and in particular, can incorporate our market-based option distributions for individual variables 
representing securities with active options markets. A slightly different approach is to use cor- 
relations of the warped variables. This procedure is likely to be more accurate, but may involve 
more computational time. 

We indicate some details. As before, it is notationally convenient to use Vi as a second notation 
for a;*, favoring the latter for fixed values and the former as a variable. Let C be the historical 
correlation matrix of the log variables (In v\ , . . . , In v n ), whose entries are the cross-correlations 

E(lnvi\nvj) — E(lnVi)E(\nvj) 

Then all diagonal terms pa are equal to 1, and C is a positive semi-definite covariance matrix, 
which we may here assume to be nonsingular (positive definite). If we use instead correlations 
of warped variables, we have simply 

Pij = E (yi2/i)- 

Let us define Fc(y\j . . . , y n ) as the cdf of a multivariate Gaussian random variable with mean 
zero and covariance matrix C. Thus, Fc(b\ , . • • , b n ) is the probability that each variable yi is 
at most some value 6*. There are more elaborate versions, such as Fc(a\ , . . . , a n ; b\ , . . . , 6 n ), 
giving the probability that each yi satisfies a$ < yi < 6*. In the single- variable case these latter 
functions are obtained from the simple cdf by a single subtraction, involving two terms, but the 
corresponding bivariate case involves four terms, and in n dimensions there would be 2 n terms. 
However, each of these more elaborate cdf's can be directly computed as an integral, just like 
the simple cdf. Since the more elaborate cdf's are needed for Monte Carlo calculations, possibly 
in high dimensions, it is best to think of them as being computed directly. 

We then define the multivariate cdf's 

F(xi , . . . , x n ) = F c (yi (x\ ), . . . , y n (xn))> 

and 

F(a x , . . . , a n \ hi , . . . , b n ) = F c {yx (ai), . . . , y n (a n ); yi (6 X ), . . . , y n {b n )) 

where the yi(xi) are the known warping functions for the individual variables. We find it conve- 
nient, with some abuse of language, to speak of F(x\ , . . . , x n ) as "the cdf" , even though we have 
all of the above functions in mind, and to use F(x\ , . . . , x n ) as a proxy for the whole distribution 
(which it does, theoretically, determine). This multivariate cdf then has the following properties: 



• Since the marginals of Fc{z\ , . . . , z n ) are Gaussian with mean 0 and variance 1, the mar- 
ginals of F(x\ , . . . , x n ) are equal to N(yi(xi)) = F(xi)] i.e., they are correct according to 
each single-variable model. 

• If the logvariables (Invi , . . . , In v n ) are actually jointly Gaussian, then the multivariate cdf 
F(x\ , . . . , x n ) is correct. 



In summary, the true joint distribution is approximated by a jointly lognormal distribution 
using historical correlations, combined with warping functions on each variable such that the 
marginal distribution of each variable is correct according to a selected single-variable model 
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(for example, according t^Wir single-variable model for optionable securities, or according to 
the lognormal model using historical volatility). The single variables may actually be portfolios, 
with a default distribution for the portfolio return being lognormal, based on historical volatility. 
This multivariate theory generalizes both our single-variable theory and standard multivariate 
(log)Gaussian models. It again allows for market input through option prices, to the extent 
that components have an active option market, but does not exclude nonoptionable securities, 
and also allows portfolios as single variables. In this way BARRA (or functionally equivalent) 
factors are also allowed because of their interpretation as portfolios of long and short positions. 



4 Applications to portfolios 



Given the multivariate cdf F(x\ , . . . , x n ) — Fc(y\ (xi), . . . , 2/ n (x n )), we can answer many typical 
questions. We first give an overview, and then take up some of the applications in more detail. 

As one example, suppose that we want to find the cdf of a portfolio variable 

X = h\X\ + fl2X2 + • ■ • + h n x n , 

where the hi are arbitrary coefficients. A simple Monte Carlo method, probably not the fastest, 
is to draw random samples from the jointly Gaussian distribution with cdf Fc(y\ , ... , y n ), trans- 
form each yi via the inverse mapping function Xi(yi), and then compute the resulting output 
sample 

x = /iixi(yi) + . . . + h n x n (y n ). 

After enough samples, we will have an approximation to the cdf of x. More precisely, the 
probability that a < x < b is, approximately, the average number of samples y\ , . . . , y n with 
a < h\xi(yi) + . . . + h n x n {y n ) < 6, and this approximation becomes exact in the limit for large 
sample sizes. This works for real portfolios, or for portfolios constructed from a number of assets 
and a residual variable, as might arise from a regression. Usually the regression is done in the 
log domain, which we discuss below. Note that the Monte Carlo method just described works 
perfectly well if the expression for x above is replaced by any function f(x\, . . . , x n ) of the x^ 
possibly quite nonlinear. 



4.1 Log domain portfolios 

In this subsection, we point out how our methods fit with another paradigm in common use in 
the financial community, and set up some further notation. It is common to work in the return 
domain, or equivalently, with logarithms; i.e., 

In x = /?i In x\ + • • * fi n In x„. 

Ignoring any possible identification of these variables with those in the previous section, the same 
discussion and Monte Carlo method as above applies, if we regard x as a nonlinear portfolio 
x = f{x\ , . . , A x n ) =^exp(/?i lnxj^ + • • • /J n \nx n ). If the sum B of the /Vs is 1, such an x may be 

written x = h\X\ + h<ix<i -1- h h n x n where hi = Pix/xi. Even if B is not 1, incremental changes 

("returns") dlnx computed from this equation for x are consistent with the above expression for 
lnx. It is common in the financial community to think of hi as approximately a constant hi, so 
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that for short periods, wn^^the x^s do not change too much, this equation for x is comparable 
to the portfolio equation in the previous subsection. 2 

For an asset x not given explicitly in terms of the terms of the x*, we obtain a similar expression 
via linear regression: 

lnx = Pq lnx 0 + Pi lnx\ + • • • + p n lnx n , 

The fii for i 7^ 0 axe correlation coefficients chosen to minimize the variance of the residual in 
historical data (perhaps subject to constraints, such as & > 0 and Yli=\ & = For example, 
x might be an security without a quoted option, and the x* for i / 0 could be taken as assets for 
which we individually know the probability distributions, in addition to the required correlation 
coefficients for x. We have written the residual term as po lnxo (usually thinking of Pq = 1 and 
the residual as normally distributed). 3 The mean of the latter could be nonzero, giving the 
regression "alpha" — a constant term making the mean of the regression correct. Alternatively, 
we could modify the equation to allow an explicit alpha, and keep the residual mean zero. 
Another minor variation might include the addition of a dummy variable with constant return, 
to adjust the value of x up or down. In particular, this gives another way of adjusting the 
residual mean to zero. This equation gives the previous one as a special case if we allow Pq = 0. 



4.1.1 Fast fits of portfolios 

One approach, which promises to be relatively fast computationally, is the following. As in the 
development of cdf-implied volatilities in Section 2, let us assume that each logvariable In x% 
above is "Gaussian" with nonconstant variance G\{xifT. In other words, the cdf is given by 
F(xi) = N(— d<i(xi,G\(xi),T)). Our aim will be to give F(x) by a similar equation, using some 
kind of fitted curve a\ (x). We will assume that we have some class of volatility curves in mind, 
with a small number of parameters which must be determined. 

If the variables lnxo, . . . , lnx n were truly jointly Gaussian, then lnx would also be Gaussian. 
Its variance would be given by the formula 



Var(lnx) = ^PiVipijPjajT, 

where pij is the correlation between lnxj and In ay, and a\T = Var(lnxj). We therefore define 
the estimate a\ (x) of a\ (x) by the conditional expectation 

<7i(x) = E{Pi<J\{xi)pijPj<Jx{xj) | lnx = y^Pi In Xj). 

i 

The calculation of the above conditional expectation may be done with Monte Carlo methods. 
In the language of nonlinear portfolios above, we would take the function /(xi, . . . ,x n ) to be 
0 outside a thin multidimensional solid enclosing the hyperplane defined by Inv = YliPi^ nv i)- 
Inside the solid we would take f(x\ , . . . , x n ) equal to the above expression for Var(ln x), divided 
by the probability of being in the solid (also a Monte Carlo calculation). In terms of samples, we 
just take the average of Var(lnx) over all the samples that end up inside the thin solid. However, 

2 Thus hi—*- = = pidhiXi, so that for small changes dxi the change dx from the first equation is 

approximately the same as would be obtained from the second. However, this relationship requires "rebalancing" 
to remain a good approximation for longer periods. 

3 For the residual term i = 0, we can use a constant variance, or impose some generic nonconstant structure 
based on observed behavior. 
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it is not necessary to compare all values of a\ (x), but only enough to tit the parameters for the 
volatility curves we are using. 

The estimated mean of In x would be In s* - \a\ (x) 2 , with s* determined as before, or replaced 
with some risk-averse estimate, to obtain the risk-averse or "true" distribution. (It is common, 
incidentally, to use factor models such as 'these to estimate a risk-averse version of Ins* = 
^jflilns* from risk- averse values of $*.) 

Also, we mention here one useful variation: We may prefer not to view the residual term 
Po In xo as part of the model, and instead write down a joint pdf only for In xi , . . . , In x n . In this 
case we can use the double expectation 

ai(x) = E{E(^2Pm(xi)pij0j(ri(xj) | lnx = ^^lnxi)), 

id i 

where the inner expectation is with respect to the variables xi,X2, . . . ,x n , and the outer ex- 
pectation is with respect to the residual. We might take the standard deviation a(xo) of the 
residual (taking (3q = 1) as a constant, determined historically, or make an estimate based on 
some leverage model. 

Now we can estimate the cdf F(x) by 

F(z) = N(-<h(xMx),T)) 

as in the univariate case. To summarize, we use our multivariate model to determine parameters 
for a univariate model of the portfolio. After that is done, we can obtain probabilities for the 
portfolio without having to go back to the multivariate model, thus achieving a savings in 
time. We could take this one step further and think of randomly generating values of cr^XjT) 
independently of any Monte Carlo philosophy (but perhaps still throwing away values of x too 
fax out-of-the-money), and then using the values obtained to do the regression required in the 
Dumas-Fleming- Whaley approach. 



5 "What-if" questions 

The multivariate distribution lends itself to the study of many questions regarding conditional 
probabilities. For example, suppose that we want to know the effect of the increase or decrease of 
some segment of the market on a portfolio, or the increase or decrease of some macro-economic 
factor. BARRA, following earlier ideas of Ross, has viewed such macro-economic factors as 
portfolios with both long and short positions. Similarly, BARRA considers market segments 
associated to price-to-earnings ratios and other fundamental parameters, as well as to industry 
groupings, as portfolios. (See the book of Grinold-Kahn cited above.) Thus, we are led simply 
to consider the effect of one portfolio on another. 

For definiteness, let us suppose the first portfolio is x, where as above 

lnx = fa lnx 0 + fi\ lnxi + • ■ • + /? n lnx n , 

and the second portfolio is y : where 

In y= 7olnyo + 7i lnxi + ••• + 7 n lnx„. 

We take /?o = 7o = 1, and view lnxo and e = Inyo as residuals with mean 0. The latter residual 
is not assumed to be a factor in our multivariate model. Consider the following typical "what-if" 



14 



l^^^en positive constants. If we know x 



question: Let A and B Driven positive constants. If we know x '^a at time T, what is the 
probability that y > B at time T? We give two approaches to this problem, the first probably 
quicker, but possibly not as accurate, using a regression to avoid at least some Monte Carlo 
calculations. 



5.1 cc What-if": An approach involving part regression, part IVfcmte Carlo 

We have In y > In B iff In y — e > In B - a. All correlations py between In xi and In xj are assumed 
known. We may also assume that we have historical values of volatilities o\ = ^/Var(ln Xi). 
(Alternatively, we could estimate such values as expected values of implied volatilities, but it 
would not be difficult to maintain an inventory of historical values, and more in the spirit of 
this part of the calculation to do so.) Thus we can estimate the historical covariances between 
lnx and Iny — e: 

Cov(lnx, Iny - e) w ^PmPitfjvfT, 
hi 

as well as a\ n x — ^/Var(lnx), a\ n y - € = -yVar(lny — e) and the correlation 

Cov(ln;r, In y — e) 



P (e) — P\nx,\ny-e — 



^ln x&\n y—e 



This gives a standard regression for the variable In y — e expressed in standard deviations from 
its mean, in terms of a similarly standardized expression for lnx. Note that e has mean 0 by 
construction. Put &i (5*, x, a) = ^]yf T ^ - Thus — (s* , x, a\ n x ) measures standardized In x 
using historical volatility, and — d%{x) = — d2{s%^x^a\{x)) measures "standardized" (warped) \nx 
using the cdf-implied volatility curve a\ (x), as discussed in the previous section. Here s* denotes 
our best estimate for the value of x at time T. 

Let <J\(y,e) denote the volatility curve associated with Iny — e, which may be estimated as 
in the previous section (or computed from estimates of a\ (y) and the standard deviation of 
the residual, if we are willing to view the residual as uncorrelated with Iny — e, as is guaran- 
teed in unconstrained regression). Put d2(y,e) = d2(s* : ye~ € ,ai(y,e)), so that — ^(y,^) is a 
"standardized" measure of Iny — e. Then the standard regression appropriate to our model is 

-<h(v>e) =P(^)(- d 2(x)). 



There is a residual associated with this regression, which we have not written down. It is pre- 
sumably normal, and its variance may be computed. For notational reasons we will just imagine 
it has been incorporated into the original e. As is apparent from the form of the expressions 
in the display, an alternative to the above regression is to do it with the warped correlation 
coefficients suggested in the previous section. If, in addition, it was appropriate to view the 
original portfolios as linear combinations of warped variables (our standard normal marginals), 
the regression above could be done without any recourse to Monte Carlo calculations. Simi- 
lar remarks would apply if we used constant historical volatility functions throughout, though 
presumably the latter procedure would lose accuracy. 

In any case, we can now answer our "what-if " question as a simple expectation in the univariate 
normal distribution of the (adjusted) residual e. Abbreviate ^(s*, A } a\ {A)) to di{A) and 
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a2\B, e)- Assume p (e) > 0 (the natural case of i 



d.2 (s*, Be £ ,&i (B, e)) to d^B, e). Assume p (e) > 0 (the natural case of a positive correlation). 
Then we have 

My>B\x>A} = E(Pi(-d 2 (y,e)>-d2(B,e)\-d 2 {x)>d2(A)) 

= E(Pr(-d 2 (x) > pie)' 1 (-MB, e)) \ - da(x) > -d 2 (A))) 
= E (min{l, Nipie)- 1 (d 2 (B, s))/N(d 2 {A))}). 

The first equation follows just because — d2(y,e) is monotonically increasing as a function of y; 
that is, the condition that y > B is completely equivalent to the condition — d 2 (y, e) > —d 2 (B ) e). 
Similar remarks hold for the condition x > A, while the expression Pr{y > B \ x > A } just means 
the probability that the condition y > B holds when it is known that x > A. The second equa- 
tion is then derived with the displayed expression above for —d 2 (y,e). (If p(e) is negative, the 
inequality involving its inverse reverses.) This inner expectation is then calculated in the normal 
distribution. For values of e for which — d 2 (A)) is as large as p(e)~ l (— d2(B } e)), the expecta- 
tion is a certainty, and yields the value 1. When — ck^A)) is smaller than p(e)~ l (— d<i(B, e)), 
its cumulative normal distribution value AT(— d 2 (^)) is smaller than N(p(e)~ 1 (—d2(B,e)), and 
the probability 1 — N(— ^(A)) = N(d2(A)) that the standard normal variable z = —d2(x) 
is at least — ^(A) is smaller than the corresponding probability 1 — N(p(e)~ l (— ^(5, e)) = 
N(p(e)- 1 (d 2 (£, e)) that z be at least p(e)" 1 (-d 2 (£, e). The ratio N{p(e)- 1 (d 2 (B, e))/N(d 2 (A)), 
which is the desired inner expectation, is thus smaller than 1, as is appropriate for a probabil- 
ity, conditional or not. If p(e) is negative, similar reasoning leads instead to the expression 
E(max{0, (N(d 2 (A)) - N(p(e)- 1 {d 2 (B,e)))/N(d 2 (A))}) for the desired conditional probability. 
Although the final answer in either case is an expectation (over e), it is essentially an integral 
that could be computed quickly with power series. (A very simple and accurate power-series 
expansion of N(z) is given on p. 252 of the book by Hull cited above.) Using that, one could de- 
termine by iterative methods what value of e makes, say, the ratio N(p(e)~ 1 (d 2 (£, e)) /iV(<2 2 (A )) 
equal to 1, and then integrate the ratio against the standard normal pdf from — oo to the deter- 
mined value of £,in the p(e) > 0 case. Similar remarks apply if p(e) < 0. (Note that, if p(e) = 0, 
the variables lnx and lny are uncorrelated, and the conditional probability Pr{y > B\ x > A} 
is the same as the unconditional probability Pr{y > S}.) 

All of the latter calculations can be done very fast. Of course, we have already used some 
Monte Carlo calculations to get this far, unless we are in the simplified context of constant 
volatility functions. 



5.2 "WhaMf": The full IVbnte 



It is easy to say how we would compute an answer to the same "what-if" question, using our 
full joint probability distribution. We simply write 

Pr{y > B\x > A} = E(Pr(-d 2 (y,e) > -d 2 (B,e) | - d 2 {x) > -d 2 (A)) 

and interpret lnx in -d 2 (x), and lny — e in —d 2 (y,e) in terms of their expansions in lnxo, 
lnxi, . . . , lnx n . To compute, say the inner expectation by a Monte Carlo calculation, we would 
generate a large number of random samples of multivariate standard normal vectors z with 
covariance matrix C. We then take the average, over the samples z which happen to satisfy 
z > —d 2 (A ), of the function which is 1 when -d 2 (y, e) > -d 2 (B, s) and 0 otherwise. We have not 
experimented to see whether this method yields better answers than the regression procedure 
above. Nevertheless, it illustrates how we could approach more sophisticated "what-if" questions 
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that could not be easily totted by regressions. For example, suppose we believe that factor w 
will remain in a range C < w < jD, and ask the same question about y y subject to the same 
condition on x. This is hard to formulate in terms of regression, and is simply not possible in 
terms of single-factor regression. However, it is easy to answer with the full distribution: 

Pr{y>B \x> A,C<w < D} = 

E(Pr(-d 2 (y,£) > -d 2 (B,e) \ -d 2 (x) > d 2 (A)), -d 2 (C) < -d 2 (w) < -d 2 (D)). 

Finally, we may not want to work in the log domain, which, if we started with a fixed portfolio 
x = h\X\ + h 2 x 2 + * • • + hnX n would force us into an approximation, as noted. But, working with 
the full distribution, we can phrase a condition x > A as hix\(yi) + . . . + h n x n (y n ) > A, in the 
language of the first section where the vector of y's plays the role of our vector z here. Monte 
Carlo calculations can now proceed as before, using log domain expressions or not for the other 
conditions. 



6 "You gotta believe" questions 

In the previous section we were focusing on an investor thinking about the value of his or her 
portfolio y in response to the change in a factor x. Conversely, an investor might want to know 
what the investment world looks like if a given stock or index y goes to a certain level B at time 
T. What is the expected value A at time T of another portfolio or simply of one of the factors 
2j? Our main plan is, upon input by the user that y is going to level B, to list several assets Xi 
or factors/indices x most highly correlated with y and their expected values with y at 5. 

It would also be possible to display a confidence interval for each selected asset or factor, 
and have other information about its new projected probability distribution readily available. 
We could also offer comparisons with the old projected probability distribution of x, where no 
assumptions on y is made. Finally, in some cases, where it was possible to explain much of 
the variance of y with just a few X{ (appearing in the regression of y), we could list percentage 
increases/decreases of a portfolio of these x% required to make B the expected value of y, based 
solely on its dependence on this portfolio. (For example, the coefficients in the portfolio could 
come from the regression of y with respect to all the x<, or some new regression might be done, 
perhaps allowing user-defined constraints). It should be mentioned that medians or modes 
are alternatives to expected values (means) here and above; in any case users will need to be 
educated about the fact that the median and mode differ systematically from the mean in near 
lognormal distributions. 

The main problem might be viewed as understanding the probability distribution of x, given 
that y > B at a given time T, with x and y as in the previous section. This can be approached 
by the methods of the previous sections, by reversing the roles of the variables. 

There is, however, a simpler question that can be treated in an especially quick way. Consider 
the problem of determining the mean of x conditioned on the equality y = B at time T. The 
idea is to use simple regression methods, but interpret answers as measured in terms of our 
variable volatilities. In our previous notation, we have a regression 

-d 2 (x) = p- (-d 2 (y,e)) + v 

where p (which we called p(e) in the previous section) is the historically determined correlation 
between In x and the random variable In y — e. Note that the roles of dependent and independent 
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variable are reversed. Therms also a residual v, which has mean 0 here, and plays no role (gets 
averaged away). Thus, the desired conditional expected value A of re is obtained from 



-d 2 (A) = E(~d 2 (x) \y=B) =p-E(-MslBe- e t ai( yi e))) = p • (-d 2 (s* y , B&fae)). 

Recall that d\{y,e) is an estimate, obtained by Monte Carlo methods, of the implied volatility 
a\ associated to the random variable Iny — s. For faster but less accurate calculations it can be 
estimated historically as YLi,jPi a iPijPj a j w ^h eac h °f the cr's, jfl's, and p's here given histori- 
cally. (See the previous section for notation.) Similarly, for fast calculations, —di{x) could use 
historical volatility, though we expect it to be given more accurately, or rather, more accurately 
according to the market view, as —d 2 {x) = —d2(s*,x,a\(x)), using the implied volatility func- 
tion estimate d\(x). If x = xi is a single asset or index in our model, then d\{x) = d\(xi) does 
not require a Monte Carlo estimate, but is presumably already available. 

To summarize, the conditional expected values required to answer "you gotta believe" questions 
are easily obtained by regression methods. The accuracy of such answers is enhanced, or at 
least shaped more to reflect market input, when all logvariables are measured in "standard 
deviations," interpreted as our variable volatilities. 

7 Portfolios containing option securities 

We conclude this document by briefly pointing out that our methods, when using full Monte 
Carlo calculations, easily apply to portfolios containing option securities. The well-known idea 
is to think of an option as as a kind of nonlinear portfolio — a quadratic one, to be more 
precise. Thus, an option on a single underlying security with underlying price x\ has a price 
approximately x = c+A(a;i — si) + (l/2)r(xi — s\) 2 for x\ near s\ , where the option was evaluated 
to a known value c. Here A and T are well-known parameters in the options markets, giving 
the first and second derivatives of the option price at s\ with respect to the underlying security 
price x\ . Perhaps the most characteristic feature of options is that they have nonzero Y — their 
proportion of increase or decrease with respect to the underlying security price changes as the 
security price changes. Explicit formulas in terms of other standard parameters are available, 
say, in the Black-Scholes theory for both A and T (see the Hull book cited above). Such formulas 
could be obtained by differentiation directly in other theories or when using empirically-fitted 
curves. In any case, once we have such an explicit approximation to x, its probability distribution 
is easily given by the Monte Carlo methods of Subsection 3.1 above. The same method applies 
as well to portfolios containing several options and other securities. 
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